%%
% int OUWdt
% index J for X
betaK=25; 
% time T
T = 1/2;
% time mesh 
num=50; tspan = linspace(0,T,num+1);
%sigma
sig=sigmahandle; a=1;b=1;m=1;sig0=1;
s  =wihandle;
IIpart=index2(betaK);
%options = odeset('AbsTol',1e-12,'RelTol',1e-6);
options = [];
row=1+2*betaK+(betaK-1)*betaK/2; X0=zeros(row,1); 
% solver of X
disp('solving X3')
tic
[time,X3] = ode45(@randouwdt,tspan,X0,options,T,a,b,m,sig0,betaK,sig,s,IIpart);
toc
%varitest(time,X3.*X3,T,'uncorrelated','oudt',num)

%%
% WCE price
Npath=10000;
xi=randn(Npath,betaK);
xisquare = 1/sqrt(2)*(xi.^2-1);
xixi=etaietaj(xi);

price=zeros(51,1);
X3=a*X3';


for iter=1:51
    asset=zeros(Npath,1);
    for ii=1:betaK
        asset=asset+X3(ii+1,iter)*xi(:,ii);
    end
    for ii=1:betaK
        asset=asset+X3(ii+betaK+1,iter)*xisquare(:,ii);
    end
    for ii=1:betaK*(betaK-1)/2
        asset=asset+X3(ii+2*betaK+1,iter)*xixi(:,ii);
    end
    asset=asset+X3(1,iter);
    price(iter)=sum(abs(asset))/Npath;
end

disp('MC price')
tic
Npath=30000;
randseed =randn(5000,Npath);
sampleou=oupath(randseed,1,1,1,1);
sampleW =oupath(randseed,0,1,1,0);
%randseed =randn(5000,Npath); % because of uncorrelated
sampleX =sum(abs(midpoint(sampleou.*sampleW)),2);
mtclprice=sampleX/Npath;
toc

hold on
plot(time,price)
plot(time,mtclprice,'r')